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Technical  Report  Summary 


Malin’s  (1980)  first-order  single  scattering  theory  has  been  extended  to  study  the 
scattering  of  surface  waves  as  well  as  body  waves  by  distributed  point  scatterers  in  a 
layered  medium.  The  scattered  waveform  itself  is  generated  and  examined  instead  of 
its  energy  envelope.  The  theory  used  allows  1)  mode  conversion  2)  wave  type 
conversion  3)  finite  scatterer  distribution,  and  4)  the  effect  of  attenuation  from 
scattering  as  well  as  intrinsic  absorption.  The  cases  studied  are  for  elastic  or  slightly 
attenuative  media  with  any  kind  of  source  and  receiver  at  any  place  in  the  layered 
structure.  This  direct  calculation  of  coda  waves  provides  us  an  immediate  description 
of  the  relation  of  coda  and  scattering.  The  objectives  are  to  find  1)  the  effect  of  layer¬ 
ing  on  scattering,  2)  the  effect  of  scatterer  distribution  on  recorded  vertical  and  hor¬ 
izontal  motion,  3)  the  relation  of  scattering  Q  to  intrinsic  Q,  4)  the  scattering  behavior 
of  surface  and  body  waves,  and  5)  the  superposition  of  scattering  waves  to  form  the 
coda.  The  generation  of  body  waves  by  ’locked  mode’  approximation,  which  makes 
the  body  wave  a  subset  of  the  'surface  wave’,  is  extensively  studied.  The  prelimilary 
results  explain  some  observed  coda  behaviors  surprisingly  well.  We  find  a  larger 
geometrical  spreading  for  near  scatterers,  which  is  caused  by  mode  conversion  or 
wavetype  conversion  because  of  the  wide  angle  scattering.  This  makes  the  spreading 
correction  higher  for  early  part  of  coda  which  may  compensate  the  low  Q  observed  in 
early  coda  of  regional  earthquakes.  This  study  is  of  practical  value  as  an  effort  to 
understand  the  complicated  coda  phases. 
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Syn thesis  of  Coda  Waves  in  the  Layered  Medium 
Chien-Ying  Wang  and  Robert  B.  Herrmann 


Introduction 

Because  of  recent  interest  on  attempting  to  extract  useful  information  from  the 
coda,  a  theoretical  study  of  coda  waves  is  urgently  required,  especially  since  the  theory 
lags  behind  observations.  Most  of  current  coda  analyses  are  based  on  Aki’s  (1969) 
backscattering  model,  which  was  extended  by  Herrmann  (1980)  to  determine  coda  Q. 
This  model  treats  the  coda  as  a  smoothly  decaying  envelope  which  is  formed  by  many 
randomly  scattered  waves.  Figure  la  is  a  real  seismogram  from  the  ECTN  network. 
The  envelope  has  an  exponential  decaying  form.  Taking  a  closer  look,  however,  the 
decay  is  not  so  uniform.  Figure  lb  shows  the  rectified  signal.  The  seismogram  looks 
as  if  it  is  composed  of  many  small  ’pulses’  or’beats.  This  is  even  more  obvious  after 
bandpass  filtering  at  1  Hz  (Figure  lc).  It  is  interesting  to  see  if  such  pulses  can  be 
analyzed  separately. 

Since  the  scattered  wave  is  the  result  of  an  averaging  process  for  waves  reacting 
with  randomly  distributed  scatterers,  most  theories  take  a  stochastic  approach  and 
describe  the  total  wavefield  by  several  statistical  quantities.  Mean  power  spectral  den¬ 
sity  is  frequently  used,  for  instance  (Aki,  1969;  Sato,  1977).  Using  such  an  approach, 
the  analysis  of  data  provides  a  gross  description  of  the  distributed  inhomogeneities  and 
the  scattering  process,  but  some  ambiguities  arise.  The  problem  of  distinguish 
between  scattering  Q  and  intrinsic  Q  is  a  typical  example  (Aki,  1982;  Dainty,  1984). 
For  our  study,  we  choose  a  deterministic  approach  to  investigate  the  effect  of  scatter¬ 
ing  on  the  generation  of  coda.  We  will  directly  examine  each  individual  scattered 
pulse  generated  in  the  layered  medium. 

Sato  (1984)  proposed  a  scattering  model  which  describes  the  body  wave  scattering 
due  to  a  plane  wave  incident  on  the  scatterers.  A  similar  model  is  also  found  in  Wu 


and  Aki  (1985).  This  model  considers  the  wave  type  conversion,  scattering  pattern, 
scatterer  size,  etc.,  but  only  for  body  waves  in  a  homogeneous  whole  space.  From 
another  view  point,  Malin  (1980)  constructed  a  scattering  model  using  surface  wave 
normal  mode  theory.  An  important  feature  of  this  model  is  that  it  is  able  to  handle 
multi-layered  medium.  The  rays  scattered  by  the  inhomogeneities  in  a  layered  struc¬ 
ture,  such  as  the  earth’s  crust,  are  numerous.  It  is  only  by  using  mode  theory  that 
scattering  in  such  a  medium  can  be  treated. 

In  this  study,  we  extend  Malin’s  first-order  single  scattering  model  to  create  the 
scattering  ’pulse’  from  scatterers  situated  in  the  layered  structure.  A  simple  review  of 
this  method  is  provided  with  the  derivation  largely  simplified.  In  order  to  include  the 
body  wave  scattering,  we  also  discuss  the  ’locked  mode’  approximation  (Harvey, 
1981),  which  makes  the  body  wave  a  subset  of  the  surface-wave  contribution.  Each 
individual  scattered  wave  is  examined  in  detail  rather  than  its  statistical  average.  The 
study  is  directed  toward  understanding  the  nature  of  the  scattered  waves  under 
different  conditions.  This  enables  us  to  analyze  the  behavior  of  each  scatterer,  and 
hence  makes  it  easier  to  isolate  important  factors  which  may  be  lost  in  an  averaging 
process. 

First-Order  Single  Scattering  Theory 

The  model  considered  is  a  layered  medium  (Figure  2]  with  inhomogeneities  dis¬ 
tributed  at  some  restricted  regions  in  it.  The  inhomogeneous  region  is  assumed  to  be 
small  compared  to  the  wavelength.  We  call  it  a  scatterer  parcel.  If  the  parcel  is  not 
small  enough  to  be  thought  of  as  a  point,  we  can  divide  it  into  several  small  parts  and 
perform  an  integration  to  find  the  total  field.  For  this  study  we  just  consider  the  point 
scatterer.  This  point  scatterer  generates  secondary  waves  when  acted  upon  by  the 
incident  wave.  The  scattered  wave  has  the  form  of  pulse  with  limited  duration  (Figure 
3).  The  coda  is  thought  to  be  composed  of  numerous  such  scattering  pulses. 
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To  satisfy  the  Born  approximation  (weak,  single  scattering),  the  inhomogeneities 
are  assumed  to  be  small  perturbations  of  the  background  material: 

P  =Po  +  Sp 
Cyu  =C0  jju  +  5C|^j 

p  is  the  density  and  is  the  elastic  constant.  For  the  isotropic  elastic  medium  con¬ 
sidered  here,  C^’s  recuce  to  the  Lame’s  constant  X  or  the  rigidity  p  .  Since  we  just 
discuss  the  scattering  from  a  small  parcel,  the  mean  behavior  of  these  inhomogeneities 
is  not  assumed.  These  three  parameters  p,  X,  and  p  can  be  reduced  to  a  single  param¬ 
eter,  if  there  exist  relations  among  them.  Sato  (1984)  used  the  velocity  perturbation; 

6a  6/3  _ 

a  ""  $ 

and  adopted  Birch’s  law  relating  the  density  and  wave  velocity 


a  =aj  p  +  ao  0  =bt  p  +  b0 

(a’s  and  b’s  are  constant).  This  results  in 


=(2  +  v)  e 


The  value  of  v  is  kept  at  0.8  in  this  report.  £  will  be  the  only  parameter  needed  to 
describe  the  inhomogeneity  variation  under  these  assumptions.  Wu  (1984)  used 
different  values  for  density  perturbation  and  elastic  constant  perturbation  which  results 
in  two  different  types  of  inhomogeneities;  impedance  type  and  velocity  type,  each  of 
which  has  different  scattering  properties.  To  reduce  the  number  of  independent 
parameters,  we  just  use  one  parameter  £  to  represent  them. 

As  indicated  by  Hudson  (1977),  the  scattered  field  for  any  first-order  scattering 
theory  can  always  be  expressed  using  the  representation  theorem: 


For  the  isotropic  elastic  case,  this  is  equivalent  to 
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uj*>  =/ dvj  bp  w2  GiJ  uf  -  b\  dp?  dp?  -  8(t  dp?  ( 3jUi0  +  dp?)  j  ( 1) 

where  u°  is  the  wave  field  incident  on  the  inhomogeneous  body,  and  G ?  is  the  Green’s 
function  for  the  corresponding  boundary  conditions.  For  the  case  of  body  wave  in  a 
homogeneous  space,  G$J  has  the  form  as  given  by  equation  (4.30)  in  Aki  and  Richards 
(1980),  which  is  used  by  Sato  (1984)  and  Wu  and  Aki  (1985).  For  the  case  of  a  sur¬ 
face  wave  in  the  layered  medium,  G ?  can  be  determined  by  the  normal  mode  theory 
(Malin,  1978). 

The  surface  wave  case  is  easier  to  construct  since  the  incident  wave  and  the  scat¬ 
tered  wave  are  trapped  in  the  layers  and  propagate  two  dimensionally.  From  the  sur¬ 
face  wave  mode  theory  (Wang,  1981),  the  incident  wave  field,  which  has  the  mode 
order  m  and  wave  type  i/  (  u  —«R  for  Rayleigh  wave  and  L  for  Love  wave),  has  the 
form 


—  - 
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where 

k  =  wave  number, 

x  =  the  distance  between  the  source  and  the  scatterer, 
a  =  amplitude  factor, 
c  =  phase  velocity, 
g  =  group  velocity, 

I0  =  energy  integral  of  surface  waves. 

U„  U„  U,  =  displacement  eigenfunctions  at  the  radial, 


(2) 


vertical  and  tangential  direction, 
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The  source  term  sm  for  double-couple  and  explosive  sources  is  given  in  the  Appendix 
I.  The  Green’s  function  with  mode  n  and  wave  type  rj  for  the  wave  field  after  scatter¬ 
ing  is 

1  -  >  ’k-r+i— 

'Gi=  1  -U,  «a,e  4  'RoJ  (3) 

V2n-  "k.r 


*U  = 

UfCos*  ,  U,sin*  ,  -iU, 

scatterer 

«R  = 

Urcos0  ,  U,sin^  ,  -  iU, 

receiver 

Lu  = 

Uisintf  ,  -Ujcos*  ,  0  ] 

*  scatterer 

LR  = 

U|sin0  ,  -Ufcos^  ,  0 

receiver 

r  is  the  distance  from  the  scatterer  to  the  receiver.  The  scattering  angle  *  and  the 
receiving  angle  <f>  are  defined  in  Figure  2a.  These  angles  are  important  since  they 
control  the  recorded  wavefield.  Scattered  P-SV  waves  can  appear  on  the  receiver  hor¬ 
izontal  component  defined  as  tangential  with  respect  to  the  source. 

Substituting  these  functions  into  the  equation  (1),  the  scattered  wave  field  1,,'umn 
which  includes  the  conversion  from  mode  m  of  wave  type  v  to  mode  n  of  wave  type  rj, 
is  written  as 


uma 


=/  dV 


„  -K  'M+ty 

a™  *a.saR,e 


2 ifyj  *Tcmx  *kBr 

(6pu2  5X  **B mu-6p  "»Cmi) 


where 


""A*,.  -Un  U„  cos 9  +  Ujm  U„ 

=(  Rkm  Ru™  -  d,u,j  (  X  Un  -  d,U,B) 

""Cm.  =2  (  “k.^cos2*  +  a.u.H.s.u.,,) 

+  (  xum  +  a,ura)  ( Rk8u„  +  d,vj  cos* 


(■*) 


^Aaj,  “-Una  U#B  sin* 

"Cm,  —  RkmUfm  LkBU,B  sin  29 

-a,u>.(dIunB+  XU,™)  Sin* 


aAM  =U,m  U„  sin* 


“Cm  =  Lk»U,m  XU™  sin 20 

+  d.U#m  (d.U™  +  XU1B)  sin0 


aAmt  =U#m  U,B  cos0  * 

LLCB1B  =  Xu^n,  LkaU»„  cos20  +  3,U|m  d,U#B  cos0 
We  also  define  a  scattering  term  FmB  as 


-iC\mx+'*kBr)+if 


(5pw2  "’A ma-S\  "'Bni-Sn  ”>Cma) 
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=/av{1 
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(5) 


]} 


The  final  total  field  solution  is  the  sum  over  all  modes  and  conversions  for  all  scatter- 
ere  in  the  region.  The  scattering  term  F  corresponds  to  the  ’volume  factor’  of  Wu  and 
Aki  (1985)  if  the  Green’s  function  related  term  D  is  not  included.  For  a  point 
scatterer,  we  simply  use 


“,F_ 


■iVe 


-i(  *krf+nBr)+ii. 


with  6V  having  the  dimension  of  unit  volume. 


Scattering  Q 

To  include  the  effect  of  energy  loss  due  to  scattering,  Malin  (1978)  considered  a 
correction  which  consevers  energy  to  the  first  order  by  choosing  a  specific  correlation 
function  for  the  scattering  cross  section.  Although  we  are  interested  in  the  behavior 
of  each  individual  point  scatterer,  the  scattering  coefficient,  or  equivalently  the  cross 
section,  is  still  needed  when  defining  the  scattering  Q. 


Let  us  consider  a  small  scatterer  parcel  containing  distributed  inhomogeneities. 
For  a  single  point  scatterer,  we  can  use  a  deterministic  scattering  formula.  However, 
when  a  wave  is  incident  on  this  scatterer  parcel,  the  size  of  the  scatterers  within  it,  the 
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relative  distance  between  individual  point  scatterers,  and  the  wavelength  of  the 
incident  wave  will  determine  the  partition  of  energy  being  scattered.  From  scattering 
theory,  the  parameter  describing  this  property  is  the  ’cross  section’  (Ishimaru,  1977, 
p.10).  Within  a  unit  of  the  scattering  angle,  the  differential  cross  section  <rd  is  defined 
by 

a  _ j.  r  (  scattered  energy  flux  density) 

d  r— oo  (  incident  energy  flux  density)  ^ 

The  total  scattering  cross  section  is  obtained  by  integration  over  all  angles: 

<r,  =f  Vide 

2r 

When  scatterers  are  distributed  with  the  mean  density  n  (particles  per  unit  area  at  a 
constant  depth  for  a  two  dimensional  case),  the  scattering  coefficient  a  which  charac¬ 
terizes  the  intensity  of  the  scattered  wave  excitation  is  given  by 


a  =n  a. 

The  reciprocal  of  a  is  usually  called  the  mean  free  path.  Scattering  coefficient  a  thus 
defined  is  equivalent  to  the  turbidity  in  the  two  dimensional  case  as  compared  to  Aki 
and  Chouet  ( 1976). 

To  obtain  the  energy  flux  density,  the  wave  field  excited  in  the  layered  medium 
should  be  taken  into  account: 


<  |u  p>  =ws  J  dzp0  <uu*> 

<>  means  the  average  over  the  ensemble.  Using  equations  (2)  and  (4),  it  is  not 
difficult  to  find  the  scattering  coefficient  for  mode  pairs  mn,  and  types  vt\  by  setting  r 
in  equation  (6)  as  the  distance  between  the  scatterer  and  the  receiver: 
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Now  we  have  the  problem  of  defining  the  average  behavior  of  the  function 
<  •  Taking  a  closer  look  at  the  scattering  region  (Figure  2b),  the  scatterers  are 

distributed  over  a  plane  area  and  the  phase  term  in  F£b  (equation  5)  is  actually 

=e-ikalx+*x7)-*a('-Tj) 

_g”‘kmX-iknr  e-‘7tkmT*-  kaTr) 

=  e-‘km’t'aV 

£B  =lcB-icm  is  called  the  Bragg  vector.  Hence,  in  evaluating  FmB  from  equation  (5), 

Fa,.  =/  dV |  e  lk“*  2  [  £(?  ,z)  e*77®  ]  Dm#( z)} 

=/dz{e  ^  ^  2  ^  (kB,z)  Daa(i)  } 

where  we  have  made  a  two  dimensional  Fourier  transform  for  the  scatterer  distribu¬ 
tion  spectra  in  the  horizontal  direction  (Kennett,  1972).  The  ensemble  average  of 
Fm  is  thus 

<F2„>  ==/ dzdz'j  <f(kB,z)  |(kB,z')>  DmB(z)  Dmn(z') 

<££>  ,  the  correlation  function  of  inhomogeneity  distribution,  is  usually  given  by  a 
form  of  an  exponential  function  or  a  Gaussian  function  (Chernov,  1960).  Assuming 
that  this  function  has  a  small  side  lobe,  i.e.,  small  correlation  distance  in  the  exponen¬ 
tial  form,  we  have  a  simpler  form: 

n<F2B>  =  [n/dzdz'  <|(kB,z)  f(kB,z')>  ]d2.(z)  =U.bD20( z) 

The  term  in  the  bracket,  which  describes  the  condition  of  scatterer  distribution  around 
the  point  considered,  will  be  kept  together  and  called  a  scattering  attenuation  unit  (U,J 
.  U,B  is  in  the  unit  of  km4.  If  UiB  is  maintained  as  a  single  value  in  the  calculations, 
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this  is  equivalent  to  stating  that  the  scattering  environment  is  similar  for  all  waves  pro¬ 
pagating  in  the  layered  medium.  The  expression  above  is  perhaps  too  simplified.  The 
correlation  function  must  be  specified  if  an  accurate  scattering  cross  section  is  desired. 
In  this  study,  we  just  use  the  simplest  form  examine  its  possible  effect. 

The  scattering  coefficient  thus  has  the  from: 

TOa  =  U„  |-  **©4  +  -i-  "*02  +  "*00 

o  2 
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0’s  are  given  in  the  Appendix  I.  We  have  used  the  distance  r  between  the  scatterer 
and  the  receiver  to  define  the  scattering  coefficient  (equation  6)  and  to  calculate  the 
incident  energy  flux  and  the  scattered  energy  flux  by  fixing  on  the  source-scatterer- 
receiver  geometry.  The  mode  pair  mn  reduces  its  energy  by  e~an*ir+x)  after  trveling  the 
distance  from  source  to  scatterer  to  receiver.  The  effect  thus  caused  by  the  scattering 
of  inhomgeneities  is  called  the  scattering  attenuation. 

Numerical  Experiment  I:  Mode-conversion  and  Scatterer  location 

With  the  velocity  model  chosen  (CUS  model,  see  Table  1),  we  generate  the 
eigenfunctions  of  surface  waves,  mainly  the  Lg  wave,  up  to  5  Hz.  Using  these,  the 
scattered  waves  are  created  by  specifying  the  source  (depth,  mechanism),  the  receiver 
(distance,  instrument),  the  scatterers  (location,  material  variation,  and  spattering 
intensity),  and  the  attenuation  condition  (intrinsic  Q  and  scattering  Q).  To  make  the 
problem  easier  to  follow,  we  fix  the  following,  except  as  indicated  otherwise: 

source  depth  —  10  km 
scatterer  depth  =  5  km 
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source  time  function  =  parabolic  with  base  0.4  second, 
receiver  distance  =  100  km 
waveform  received  —  ground  velocity 
in  homogeneity  variation  =  5% 

mode  conversion  =  up  to  10  neighboring  modes  (due  to  computer 
restriction  on  speed  and  storage) 
wave  type  conversion  =  yes 

An  arbitrary  source  mechanism  was  used  and  the  effect  of  its  radiation  pattern  is  aver¬ 
aged  out  by  distributing  several  scatterers  on  the  ellipse.  In  this  section,  we  will  dis¬ 
cuss  the  properties  of  the  scattered  waves  under  different  scattering  conditions  but 
ignore  the  problem  of  scattering  and  intrinsic  attenuation  at  this  moment.  The  so 
called  ’diagonal  selection’  rule  (Malin,  1980)  which  states  that  a  particular  mode  of 
surface  wave  mainly  scatters  into  the  same  mode  without  significant  conversion  will  be 
examined  first.  Sscattering  from  different  locations,  e.g.,  distances,  depths,  scattering 
angle,  etc.,  are  also  discussed. 

Malin  (1980)  used  scattering  coefficient  of  vertical  component  from  an  acoustic 
model  to  justify  the  diagonal  selection  rule.  Since  the  scattered  wave  instead  of  its 
energy  envelope  is  generated  here,  this  rule  will  be  discussed  by  examining  the  resul¬ 
tant  waveforms.  Two  cases  are  studied;  one  is  for  the  wave  type  conversion  (from  R 
to  L  and  L  to  R),  and  the  other  for  the  mode  conversion.  Figure  3  shows  the 
waveforms  of  signals  directly  from  the  source  (  surface  wave)  or  from  a  scatterer 
(scattered  waves)  after  traveling  the  same  distances  as  indicated  at  the  ends  of  seismo¬ 
grams.  Scattered  waves  have  more  high  frequencies  and  are  of  shorter  duration  than 
the  surface  waves.  The  major  difference  is  in  the  excitation  of  the  fundamental  mode. 
The  fundamental  mode  excited  by  the  source  does  not  excite  the  scatterers  as  well  at 
high  frequencies  as  do  the  higher  modes.  This  is  a  function  of  the  eigenfunction  distri¬ 
bution  with  depth. 


Figure  4  examines  the  effect  of  wave  type  conversion  on  scattering.  In  this  figure 
and  those  that  follow,  the  number  of  kilometers  indicates  the  total  scattering  distance, 
i.e.  the  distance  from  the  source  to  the  scatterer  plus  the  distance  from  the  scatterer  to 
the  receiver.  The  distance  from  the  source  to  the  receiver  is  fixed  at  100  km.  The 
angle  in  degrees  is  the  location  of  the  scatterer  on  an  ellipse  which  has  the  source  and 
the  receiver  at  its  foci  and  the  coordinate  origin  at  its  center.  Zero  degrees  indicate  the 
scatterer  is  behind  the  receiver,  and  90  degrees  indicate  that  the  scatterer  is  equidistant 
from  the  source  and  receiver.  Three  different  components:  vertical  Z,  radial  R  (with 
respect  to  source-receiver  coordinate),  and  the  tangential  T  are  also  indicated.  In  each 
pair  of  traces  of  Figure  4,  the  lower  one  includes  the  wave  type  conversion  and  the 
upper  one  does  not.  For  each  distance,  each  pair  of  traces  is  plotted  using  the  same 
scale.  From  this  comparison,  we  find  that  the  wave  type  conversion  is  not  obvious 
except  at  short  distance  for  wide  angle  scattering.  The  Z  component  seems  to  be  more 
independent  of  this  effect  than  the  other  components.  Thus  it  is  concluded  that  the 
wave  type  conversion  is  important  only  for  the  horizontal  components  at  short  dis¬ 
tances.  It  is  at  these  distances  that  wide  angle  scattering  occurs.  For  large  scattering 
distances  relative  to  the  source-receiver  distance,  backscattering  is  all  that  occurs  from 
the  the  geometry. 

Figure  5  illustrates  the  effect  of  mode  conversion.  Three  cases  are  displayed  in 
which  the  number  of  modes  permitted  to  be  converted  to  other  near  neighboring 
modes  are  1,  5  and  10.  In  addition,  the  effects  of  wavetype  conversion  are  also  shown. 
The  scattered  waveforms  are  created  by  combining  the  effect  of  seven  scatterers 
evenly  distributed  on  the  ellipse  defined  by  the  scattering  distances.  The  amplitude 
scale  of  plotting  for  every  scattering  distance  is  the  same.  As  in  Figure  4,  the  Z  com¬ 
ponent  seems  stable  and  does  not  accumulate  as  much  contribution  from  other  modes 
except  at  the  short  distance  where  mode  conversion  must  still  be  considered.  On  the 
other  hand,  the  R  and  T  components  (Figures  5b  and  5c)  are  sensitive  to  the  number 
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of  modes  allowed  to  convert  each  other  and  to  the  wave  type  conversion.  A  low 
mode,  low  frequency  motion  can  be  seen  (Figures  5b  and  5c),  which  is  caused  by  the 
conversion  of  higher  modes  to  lower  modes  during  the  scattering.  These  results  do 
not  support  the  diagonal  selection  rule,  although  the  use  of  it  alone  will  not  cause 
major  differences  in  the  coda  excitation  For  the  elastic  medium,  the  energy  exchange 
by  the  mode  conversion  is  important  especially  for  the  two  horizontal  components. 
The  mode  conversion  can  be  looked  upon  analagous  to  a  change  of  incident  angle  on 
the  layered  medium  after  scattering.  It  describes  the  effect  similar  to  the  scattering 
pattern  for  body  body  wave  illustrated  in  Aki  and  Richards  (1980,  chap  13).  Figure  6 
is  a  similar  display  as  Figure  5  but  for  an  explosive  source  with  the  scattering  distance 
at  110  km.  A  noticeable  tangential  component  is  obtained  which  is  simply  a  result  of 
scatterer  geometry.  Mode  conversion  and  wave  type  conversion  are  both  important  at 
this  scattering  distance. 

The  location  of  the  scatterer  at  different  positions  with  respect  to  the  source  and 
the  receiver  is  important  in  exciting  the  scattered  wave.  Figure  7  is  for  scatterers  at 
different  locations  on  an  ellipse  at  scattering  distances  of  110,  300,  and  500  km.  The 
amplitude  of  each  trace  is  normalized  to  the  same  distance  using  a  l/r™  (m=l.O)  rule 
which  will  be  discussed  later.  In  these  figures,  waves  from  three  different  double¬ 
couple  sources  are  summed  to  average  out  the  radiation  pattern  effect.  We  can  see 
that  the  amplitude  changes  along  different  azimuth.  An  amplitude  relation  between  R 
and  T  components  is  observed,  e.g.,  when  R  is  large  T  is  small.  It  is  worth  noting  that 
a  scatterer  close  to  the  source  or  the  receiver  is  not  necessary  to  create  large  scattering. 
From  the  variation  of  amplitudes  around  the  ellipse,  we  can  see  an  averaged  effect  of 
scattering  pattern. 

We  next  examine  the  effect  at  different  vertical  depths.  A  result  is  shown  in  Fig¬ 
ure  8.  The  displays  are  for  scatterers  at  the  distances  200  km  (Figure  8a)  and  500  km 
(Figure  8b).  The  effect  of  two  source  depth  10  km  and  25  km  corresponding  to  the 
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upper  and  lower  sets  of  traces  are  shown.  A  Q  model  10U0.5Q,  which  will  be 
explained  later,  is  also  applied.  It  is  found  that  the  scatterer  at  greater  depth  preferen¬ 
tially  excites  the  scattered  waves  with  the  higher  phase  velocities.  The  scatterer  in  the 
upper  sedimentary  layer  (0.5  km)  is  capable  of  generating  well  dispersed,  long  duration 
scattering  waves.  Because  of  dissipative  attenuation,  the  shallow  scatterers  yield 
waveforms  that  lose  their  energy  and  frequency  content  faster,  which  makes  them 
more  comparable  to  the  the  effect  of  deeper  scatterers.  At  larger  propagation  dis¬ 
tances,  the  difference  between  scattered  waves  at  different  depths  is  gradually  reduced, 
although  the  effect  of  attenuation  may  still  be  seen.  In  Figure  8b,  we  can  see  that  the 
high  frequency  components  are  attenuated  for  shallow  scatterers  where  Q  is  low  after 
traveling  long  distance.  But  their  amplitudes  and  waveforms  are  not  too  much 
different  except  for  the  low  frequency  fundamental  mode.  This  indicates  that  a 
scatterer  at  any  depth  in  the  crust  generates  a  similar  scattered  signal. 

Numerical  experiment  II:  Scattering  Q  and  Intrinsic  Q 

As  shown  in  the  previous  section,  the  scattering  coefficient  is  a  function  of  the 
correlation  between  scatterers,  the  inhomogeneity  size,  and  the  incident  wave  length. 
To  make  the  problem  easier  to  solve,  we  defined  a  scattering  attenuation  unit  (U„) 
from  which  we  obtained  the  scattering  coefficient  a  after  taking  into  account  the  angu¬ 
lar  scattering  variation.  The  amplitude  decay  is  given  by  e~or/'!  where  r  is  the  scatter¬ 
ing  distance.  In  fact,  this  is  a  definition  of  turbidity  which  is  used  to  describe  the 
scattering  energy  loss  mechanism  (Dainty,  1981).  In  our  model,  the  scattering 
attenuation  is  assumed  the  same  for  different  scatterers  which  are  distributed  ran¬ 
domly  and  uniformly  along  scattering  wave  path  passing  through  them.  We  will  exam¬ 
ine  the  effect  of  scattering  on  the  wave  energy  loss  and  compare  this  with  the  attenua¬ 
tion  due  to  anelastic  absorption.  The  anelastic  attenuation  is  calculated  by  perturbing 
the  elastic  constants  and  combining  with  the  intrinsic  Q  values. 
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Figure  9  shows  the  effect  of  different  levels  of  scattering  attenuation  from  scatter¬ 
ed  at  a  5  km  depth  and  a  110  km  scattering  distance.  Three  plots,  a,  b,  and  c 
correspond  to  10,  100,  and  1000  scattering  units  (U,»).  It  is  found  that  high  frequency 
signals  are  just  partially  absorbed  for  different  U,»  values  at  this  short  propagation  dis¬ 
tance.  The  effect  of  scattering  attenuation  seems  to  broaden  the  pulse,  i.e.,  redistri¬ 
bute  the  wave  energy  over  longer  time  range. 

Figure  10  shows  an  example  with  scattering  distance  extended  to  300  km.  Figure 
10a  is  for  the  scattered  at  5  km  depth  and  10b  for  the  scattered  at  0.5  km  depth.  The 
uppermost  trace  in  this  figure  is  the  reference  trace  which  does  not  involve  any  Q 
effect.  Other  traces  are  named  by  the  symbol  pUqQ  which  is  read  as  scattering  Q  hav¬ 
ing  p  U„  and  intrinsic  Q  having  values  from  the  model  q.  The  basic  Q  model  used  is 
Q f  =  300  for  upper  20  km  and  2000  for  the  rest  in  the  CUS  model.  2Q  indicates 
twice  as  much  attenuation  as  IQ,  e.g.,  Qtf=450.  It  is  obvious  that  the  scattering  Q 
reduces  the  high  frequency  energy  less  than  intrinsic  Q.  This  means  that  the  scatter¬ 
ing  Q  is  not  too  sensitive  to  the  scattering  environment.  If  the  anelastic  part  of  struc¬ 
ture  does  not  dissipate  all  high  frequency  signal,  the  scattering  Q  may  not  reduce  it 
even  under  strong  scattering  conditions.  Furthermore  the  scattering  Q  reduces  the 
scattered  pulse  amplitude  but  does  not  affect  its  frequency  content  much  as  can  be 
seen  from  the  traces  in  10U0Q  and  40U0Q.  In  Figure  10b,  we  also  find  that  the 
scattering  Q  seems  to  suppress  the  lower  mode  signals  and  enhance  the  high  frequency 
higher  modes.  However,  this  feature  probably  depends  on  the  velocity  model  used. 
From  this  test,  we  conclude  that  the  effect  of  scattering  Q  and  intrinsic  Q  are 
apparently  different. 

If  the  scattered  arrivals  from  different  distances  constitute  the  coda  at  different 
lapse  time,  the  frequency  content  of  these  pulses  may  indicate  the  effect  of  Q  at 
different  places  of  coda.  To  examine  this  effect,  we  generate  the  scattered  pulses  at 
two  distances:  200  km  and  500  km  under  different  Q  conditions.  These  are  shown  in 
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Figures  11  and  12.  The  Parts  a  and  b  are  for  the  scatterer  depths  at  5  km  and  15  km 
at  these  depths  the  Q;  values  are  different.  The  spectra  of  the  traces  are  also  plotted, 
shifted  by  a  factor  of  ten  for  clarity.  One  thing  revealed  from  this  test  is  that  the 
scattering  attenuation  absorbs  the  high  frequency  signal  equally  well  irrespective  of  the 
scattering  distance  or  the  scatterer  depth.  This  is  shown  by  the  difference  in  the  slopes 
between  the  first  and  the  second  spectra.  This  is  due  to  the  fact  that  a  similar  scatter¬ 
ing  environment  is  assumed  when  calculating  scattering  coefficient  a.  The  intrinsic  Q, 
on  the  other  h-nd,  attenuates  the  high  frequency  signals  rapidly  for  the  shallower 
scatterer  where  intrinsic  Q  is  low.  Figures  12a  and  12b  can  be  thought  of  as  two 
extreme  cases:  in  12a  intrinsic  Q  dominates,  and  in  12b  scattering  Q  dominates.  The 
final  Q  value  must  result  from  the  interaction  of  these  two  attenuation  mechanisms. 
If  the  data  show  a  strong  frequency  dependence  at  high  frequency  range,  it  is  possible 
this  is  the  attenuation  from  deep  structure  and  the  effect  of  scattering  Q  prevails.  At 
the  low  frequency  (<1  Hz),  it  is  difficult  to  tell  the  difference  between  these  two 
mechanisms. 

Numerical  Experiment  III:  Geometrical  Spreading 

Aki  (1969)  proposed  a  scattering  theory  which  describes  the  coda  waves  of  a  local 
event  as  being  composed  of  backscattered  waves  from  many  random  distributed  inho¬ 
mogeneities  in  the  lithosphere.  A  relation  from  his  model,  which  has  been  widely 
used  in  the  literature,  is 

A(w|t)  =-C(w)t-‘'e*“v5Qc 

where  A(w|t)  is  the  average  peak-to-peak  amplitude  around  the  frequency  w  and  the 
time  t,  C(w)  is  the  coda  shape  function  which  includes  the  scattered  wave  excitation 
term  and  a  dispersion  correction,  f1'  is  for  geometrical  spreading,  and  the  exponential 
term  describes  the  dissipative  attenuation.  The  value  of  v  is  taken  as  1.0  for  body 
waves,  0.5  for  surface  waves,  and  0.75  for  a  diffusion  model  (Aki  and  Chouet,  1976). 
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Treating  coda  of  local  earthquakes  (distance  <  100  km),  the  source  and  the  receiver 
are  set  at  the  same  place  to  simplify  the  calculation  of  backscattering  time  in  this 
model.  Sato  (1977)  proposed  a  modification  of  the  geometrical  spreading  factor  for 
near  scatterers  which  create  waves  arriving  at  the  receiver  just  after  the  direct  wave. 
This  modification  takes  into  account  the  separation  between  the  source  and  the 
receiver  and  is  used  for  regional  earthquakes  (Pulli,  1984).  Kopnichev  (1977)  also 
derived  a  similar  correction  by  considering  the  ellipticity  of  the  scatterer  surface  for 
arrivals  at  a  given  time.  In  fact,  this  correction  represents  a  geometrical  change  of  the 
scatterer  distribution  ranging  from  a  flattened  ellipsoid  to  a  shape  similar  to  a  sphere  at 
large  distances  (from  an  ellipse  to  a  circle  for  surface  waves  as  illustrated  in  Figure 
2a).  It  is  purely  a  geometrical  adjustment.  However,  as  discussed  in  preceding  sec¬ 
tions,  the  scatterers  at  the  near  distances  yield  significant  mode  conversion  and  wave 
type  conversion.  A  more  complicated  nature  of  scattered  waves  is  expected  there. 
Similarly,  because  of  the  importance  of  scattering  angle  in  controlling  the  scattering 
mechanism  (Wu,  1984),  its  effect  must  be  significant  at  short  distance  where  the 
scattering  angle  changes  abruptly.  In  this  section,  we  will  measure  this  effect  from 
simple  nuerical  studies.  The  result  indicates  that  the  so  called  ’Sato’s  correction’ 
which  gives  a  heavier  geometrical  spreading  correction  in  the  early  part  of  coda  than 
Aki’s  formula  is  still  too  small  as  the  scattered  wave  energy  changes  rapidly  between 
modes  and  wavetypes  at  near  distances. 

In  order  to  see  the  variation  of  geometrical  spreading  factor,  we  generate  many 
scattered  pulses  at  different  scattering  distances  ,  ignoring  dissipative  attenuation,  both 
of  intrinsic  and  scattering.  An  averaged  peak-to-peak  amplitude  is  then  measured  , 
with  narrow  band  filtering  if  needed,  and  plotted  in  log-log  scale  with  the  distance. 
One  of  the  results  is  shown  in  Figure  13a.  The  scattering  distance  is  spacing  by  10  km 
up  to  600  km  and  the  receiver  is  kept  at  100  km  from  the  source.  It  is  obvious  that 
there  exists  two  slopes.  A  steeper  slope  with  amplitude  d'  upping  rapidly  with  the  dis- 
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tance  is  seen  in  the  front  part.  A  similar  tendency  is  also  observed  for  the  tangential 
component  (Figure  13b).  In  these  figures,  the  unfiltered  scattered  wave  (upper  one) 
and  the  filtered  scattered  wave  all  follow  similar  trends.  This  is  due  to  a  similar  fre¬ 
quency  response  for  scattered  wave  even  traveling  over  different  distances  (see  Figure 
11).  Figure  14  shows  a  similar  amplitude  decay  for  the  west  United  States  model 
(Table  1).  Here  we  have  set  the  source-receiver  distance  at  300  km  (Figure  14a)  and 
at  100  km  (Figure  14b).  The  rapid  amplitude  decay  happens  at  the  scatterer  distances 
smaller  than  about  600  km  and  200  km  respectively.  Hence,  we  may  define  the  region 
with  scattering  distance  smaller  than  twice  the  source- receiver  distance  as  an  ’unstable 
region’  which  has  a  rapid  decrease  in  amplitude,  and  a  ’stable  region’  at  larger  dis¬ 
tances  which  has  a  spreading  factor  close  to  one.  This  distance  separation  agrees  with 
that  of  ’Sato  correction’  of  twice  the  source-receiver  distance.  Note  that  the  decay  rate 
in  unstable  region  is  greater  than  that  in  stable  region  by  0.5-1.0.  This  value  is  larger 
in  WUS  model  than  in  CUS  model  (Table  1).  Taking  this  into  account,  we  propose  a 
geometrical  spreading  correction  basing  on  an  extension  of  Sato’s  (1977)  model  to  the 
surface  wave: 


where  a  is  the  ratio  of  lapse  time  on  the  coda  to  the  Lg  wave  propagation  time  and  a  is 
between  0.5  and  1.0.  This  correction  is  plotted  in  Figure  15. 

In  propagating  in  a  layered  medium,  not  only  does  the  amplitude  of  the  scattered 
pulses  decay  with  the  distance  but  its  duration  also  increases  because  of  surface  wave 
dispersion.  Figure  18a  displays  the  scattered  pulses  at  different  distances  and  Figure 
16b  is  a  linear  regression  which  shows  the  rate  of  duration  expansion.  The  CUS  model 
exhibits  an  approximate  increase  in  duration  of  0.05  sec/km  while  WUS  model  exhi¬ 
bits  an  increase  of  0.5  sec/km.  Because  a  thicker  sedimentary  layer,  the  WUS  model 
is  said  to  be  more  dispersive.  This  dispersion  induced  pulse  duration  expansion  must 
be  considered  when  summing  all  pulses  to  form  the  coda. 
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Numerical  Experiment  IV:  Superposition  of  Scattered  Pulses 

We  have  discussed  the  properties  of  a  single  scattered  pulse.  The  next  question 
is  how  these  pulses  are  combined  to  form  the  coda.  In  this  section,  we  will  conduct 
an  interesting  experiment  using  the  parameters  determined  in  the  preceding  sections 
to  synthesize  the  coda  envelope.  From  recent  studies  of  coda  Q  of  regional  earth¬ 
quakes,  Pulli  (1984)  and  Shin  and  Herrmann  (1985)  both  observed  a  lower  Q  values 
from  the  early  part  of  coda.  Der  et  al.  (1984)  also  mentioned  an  obvious  coherence 
difference  between  early  coda  and  later  coda.  This  section  will  attempt  to  propose  an 
explanation  to  these  observations. 

From  the  study  of  scattered  pulses  just  discussed,  we  make  the  following  assump¬ 
tions;  1)  the  shape  of  scattered  pulse  can  be  approximated  by  an  exponential  function 
te"**;  2)  the  width  of  the  pulse  depends  on  its  traveling  distance  following 

Tft 

duration^, 'distance  ;  3)  the  maximum  amplitude  of  the  pulse  decays  as  e  Q  /tm,  where 
the  exponential  term  describes  the  effect  of  attenuation  including  intrinsic  or  scattering 
Q,  and  l/t“  is  the  geometrical  spreading  factor,  4)  the  pulses  arrive  at  the  receiver  in  a 
random  manner,  thus  the  number  of  pulses  in  a  time  interval  at  any  time  of  seismo¬ 
gram  is  assumed  to  be  approximately  constant.  To  take  into  account  reasonable  varia¬ 
tion,  we  also  allow  a  10% random  variation  for  the  parameters:  m  (spreading  factor),  ao 
(duration  expansion)  and  dtp  (pulse  interval).  We  will  discuss  the  pulse  at  1  Hz  with 
the  propagation  velocity  3.5  km/sec. 

With  each  pulse  defined,  the  coda  is  then  synthesized  by  summation  of  many 
pulses.  One  example  of  synthetic  coda  thus  formed  can  be  seen  from  Figure  16  where 
every  individual  pulse  and  the  coda  envelope  after  summation  are  plotted  assuming 
different  parameters.  Since  the  energy  of  random  arrivals  is  mixed,  we  use  root- 
mean-square  rule  to  form  the  coda  enevlope.  After  constructing  the  coda,  we  then 
apply  the  Aki-Chouet’s  (1978)  model  with  surface  wave  spreading  correction  ,  i.e. , 
1/t0-6,  to  calculate  the  coda  Q  value.  The  Q  values  at  two  regions,  first  100  second  and 
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the  rest,  are  calculated  separately.  Here  a  sliding  window  with  10  seconds  length 
(Pulli,  1984)  is  also  used  when  taking  points  from  the  coda  shape.  The  cases  for 
different  values  of  m,  a*  and  dtp,  which  represent  different  scattering  conditions,  are 
examined.  The  coda  Q  values  thus  obtained  are  summaried  in  Table  2  when  theoreti¬ 
cal  Q  value  is  kept  at  1000  and  in  Taole  3  when  Q  is  150. 

One  astonished  result  revealed  by  this  simple  test  is  that  the  pulse  spreading  fac¬ 
tor,  m,  is  the  most  important  factor  which  controls  the  coda  shape  decay.  The  dura¬ 
tion  factor  and  the  pulse  density,  on  the  other  hand,  only  have  limited  effect  on  the 
coda  formation  if  the  pulses  are  crowded  together  enough  to  make  a  smoothly  decay¬ 
ing  coda  shape  as  usually  observed.  The  dispersion  of  scattered  pulses  in  the  layered 
medium,  which  causes  the  duration  increase,  dose  not  affect  the  coda  as  much  as  first 
expected  (Dainty,  1984).  We  also  found  that  for  a  low  Q  region  (Table  3)  the  coda 
decays  rapidly,  thus  not  too  many  pulses  are  needed  and  its  Q  value  comes  out  fairly 
stable.  The  coda  is  said  to  be  ’saturated’  by  attenuation.  For  high  Q  region  (Table  2), 
different  pulse  parameter  will  play  an  important  role  on  the  coda  formation. 

Two  different  Q  estimates  are  also  observed  in  this  test.  The  early  part  of  coda 
always  has  a  lower  Q  value.  Since  the  pulses  in  the  front  have  high  amplitude  and 
short  duration  and  pulses  in  the  later  part  are  lower  and  longer,  the  decay  of  the  coda 
shape  composed  of  these  pulses  must  be  more  rapid  in  the  early  part.  This  effect  is 
due  to  pulse  duration  increase  induced  by  dispersion.  This  is  another  factor  affecting 
the  early  part  coda  shape,  which  makes  the  ’unstable  region’  defined  before  even 
more  unstable. 

As  a  result  of  repeating  the  above  experiment  for  different  combination  of 
parameters  m,  ao,  and  dtp,  we  obtained  a  prelimilary  understanding  of  the  constitution 
of  coda  by  the  summation  of  scattering  waves: 

1)  For  a  high  Q  region,  such  as  the  eastern  United  States,  the  manner  scattered 
pulse  decay  with  the  distance,  because  of  geometric  spreading  or  mode  conversion, 
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controls  the  value  of  Q  determined  from  coda  shape.  A  large  m  value  in  the  early  part 
of  coda  (unstable  region)  has  a  major  effect  on  Q  value.  This  is  the  case  predicted  in 
Figure  17a. 

2)  For  a  low  Q  region,  such  as  the  western  United  States,  the  pulses  which  arrive  at 
the  receiver  to  form  the  coda  are  strongly  affected  Vv  the  attenuation  and  thus  the 
coda  shape  is  a  stable  estimator  of  Q.  This  is  shown  in  Figure  17b. 

3)  The  pulses  with  longer  duration  yield  a  longer  tail  which  in  turn  leads  to  a  higher 
Q  estimate,  but  this  effect  is  not  so  outstanding  as  expected. 

4)  If  the  pulses  do  mix  to  form  the  coda,  the  density  of  pulse  arrivals  is  not  very 
important  with  respect  to  the  final  coda  Q  values.  This  is  due  to  the  fact  that  the 
pulses  are  superposed  by  means  of  the  energy  which  reduces  the  difference  in  the 
number  of  pulses. 

5)  Estimation  of  coda  Q  using  only  earlier  parts  of  the  coda  requires  knowledge  of 
the  rate  of  correct  geometrical  spreading  term  there. 

Numerical  Experiment  V:  Body  Wave  Scattering 

The  theory  of  body  wave  scattering  in  the  earth  has  been  widely  discussed  from 
the  points  of  view  of  scattering  from  a  discrete  volmue  (  Knopoff,  1959;  Wu,  1984)  or 
superposition  of  scattered  body  waves  (Sato,  1984).  These  theories,  however,  use  a 
plane  wave  reacting  with  the  scatterer  in  a  whole  space  or  at  most  a  halfspace.  In  this 
section,  we  will  discuss  body  wave  scattering  in  the  layered  medium.  The  layered 
earth  structure  is  an  important  factor  which  provides  a  circumstance  for  waves  to  react 
with  inhomogeneities  by  trapping  waves  in  the  layers.  We  will  simply  extend  the 
theory  developed  for  surface  waves,  i.e. ,  Lg  waves,  to  body  waves  by  using  ’locked 
mode’  approximation  (Harvey,  1981).  This  is  done  by  adding  a  rigid  cap  layer  at  great 
depth  in  the  original  model  to  trap  the  body  waves.  In  doing  this,  the  body  waves 
become  a  subset  of  the  surface-wave  filed.  Hence,  we  can  easily  extend  the  first-order 
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scattering  theory  to  cover  this  part  of  scattering  signals. 

If  the  technical  difficulties  mentioned  in  the  Appendix  II  can  be  overcome,  the 
body  wave  scattering  waveform  will  be  generated  as  shown  in  Figure  18.  Three  com¬ 
ponents  of  ground  motion  with  scatterer  at  110  km  distance  and  5  km  depth  in  the 
CUS  model  are  displayed.  The  scattered  waveform  for  surface  wave  is  also  plotted  for 
comparison.  As  expected,  the  scattered  waveform  still  resembles  the  shape  of  direct 
waves.  Body  wave  scattering  is  strong  around  the  P  and  S  arrivals.  This  is  probably 
caused  by  the  P  to  P  and  P  to  S  scattering  after  the  P  wave  and  S  to  S  scattering  after 
the  S  wave.  The  scattering  from  S  to  P  seems  weak.  We  can  also  see  some  body-wave 
and  surface-wave  scattering  interaction  in  front  of  or  after  the  main  pulses.  Since  the 
amplitude  of  scattered  body  waves  is  much  smaller  than  scattered  surface  waves,  we 
extract  the  body  wave  scattering  part  by  windowing  and  display  them  at  a  larger  scale 
as  shown  in  Figure  19.  Results  at  three  scattering  distances  are  presented.  The  body 
wave  scattering  is  complicated  for  this  layered  structure.  The  shear  scattering  is  strong 
on  the  two  horizontal  components.  Because  of  the  compicated  scattered  waveform,  it 
is  easier  for  the  generation  of  body  wave  coda  especially  at  high  frequencies  after 
superposition  of  such  pulses. 

Next  we  check  the  body  wave  scattering  under  different  attenuation  condition.  In 
Figure  20a,  we  plot  the  scattered  pulses  including  body  waves,  surface  waves  and  their 
spectra.  The  first  trace  is  the  pulses  from  full  wave  scattering,  the  second  is  a 
representation  of  the  body  wave  part  from  the  first  trace,  and  the  third  trace  is  the 
surface  wave  scattered  pulse  which  is  generated  by  using  crustal  surface-wave  modes 
only.  Their  spectra  are  plotted  with  the  third  (surface  wave)  two  amplitude  order 
lower.  After  applying  Q  models  due  to  scattering  attenuation  or  intrinsic  attenuation, 
we  obtain  the  attenuated  scattering  pulses  as  shown  in  Figures  20b  and  20c.  It  is 
interesting  to  see  that  the  body  wave  scattering  suffers  much  less  attenuation  than  sur¬ 
face  wave  scattering  especially  in  high  frequency  range.  This  can  be  seen  more 


apparently  in  Figure  21  where  spectra  under  different  Q  conditions  are  plotted  for  sur¬ 
face  and  body  wave  scattering  respectively.  Figure  21a  is  for  scattering  distance  at  250 
km  and  21b  at  500  km.  In  these  figures  the  high  frequency  spectrum  slopes  do  not 
show  as  much  variation  for  body-wave  scattering  as  for  surface-wave  scattering.  The 
body-wave  scattering  contains  higher  frequencies  initially  and  its  energy  is  not  much 
attenuated  by  either  intrinsic  or  scattering  attenuation.  This  might  make  the  body 
wave  scattering  dominant  at  high  frequencies  as  observed  by  Shin  and  Herrmann 
(1985). 

Conclusion 

After  examing  the  behavior  of  scattered  pulses  generated  in  a  layered  medium, 
we  have  the  following  conclusions: 

1)  The  ’diagonal  selection’  rule  which  predicts  no  mode  conversion  is  not  proper  for 
scattering  in  the  layered  elastic  medium.  Mode  conversion  must  be  considered  espe¬ 
cially  for  the  horizontal  components  or  at  short  distances. 

2)  The  wave  type  conversion  between  Rayleigh  and  Love  waves,  on  the  other  hand, 
can  be  ignored  except  at  short  distances. 

3)  The  scatterers  near  the  surface  can  generate  well  dispersed  scattered  waves,  but 
they  also  lose  energy  faster  because  of  dissipative  attenuation.  The  combined  effct  of 
attenuation  and  layering  makes  the  scatterer  depth  dependence  less  significant. 

4)  The  variation  of  scattering  Q  does  not  affect  the  value  of  apparent  Q  as  much  as 
the  intrinsic  Q.  For  the  same  scattering  environment,  the  attenuation  of  high  frequen¬ 
cies  from  scattering  Q  does  not  seem  to  depend  much  on  the  scattering  location 

5)  It  is  difficult  to  make  a  distinction  between  the  effect  of  scattering  Q  and  intrinsic 
Q  on  the  final  Q  value. 

6)  Larger  geometrical  spreading  decay  is  found  in  the  early  part  of  coda.  It  falls  off 
by  a  factor  of  tr0,6  to  t‘l  faster  than  the  later  part.  The  geometrical  correction  in  the 
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early  unstable  region  of  coda  is  suggested  to  be  larger  than  the  simple  Sato’s  correc¬ 
tion.  This  is  caused  mainly  by  rapid  mode  or  wave  type  conversion  from  scatterers  at 
near  distances. 

7)  If  the  scattering  pulses  which  constitute  the  coda  are  dense  enough,  which  is 
required  to  make  a  smooth  decaying  coda  envelope,  the  geometrical  correction  factor 
will  control  the  value  of  coda  Q. 

8)  The  body-wave  scattering  internally  yields  higher  frequency  arrival  than  the 
surface-wave  scattering.  It  also  suffers  less  attenuative  absorption  than  the  slower  sur¬ 
face  waves. 

9)  The  effect  of  layering  on  the  coda  formation  is  very  important.  Because  of  wave 
trapping  in  the  layers,  even  a  weak  scatterer  distribution  may  create  oberverable  coda. 
This  layering  is  the  same  reason  that  makes  the  Lg  phase  the  most  dominant  phase  in 
regional  earthquakes  seismograms. 
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Appendix  I 


I.  Source  Function: 


For  a  double  couple  source,  the  source  mechansim  can  be  determined  by  the  dip 
d  ,  slip  s,  and  strike  <f>  in  the  following  forms: 

"S  *  *k  U^,  +  ( Rk  Ur)  R„  +  i(  Rk  U,+^L)  R* 


=  Lk  U,R^  -  i-^-Ri, 

dz 


where 


R„  =-sind  coss  sin2 <f>  -  -£-sin2d  sins  cos2 4> 

it 

Raa  =sins  sin2d 

Ra,  =  cos2d  sins  sin^  -  cosd  coss  cos0 

R,'a  =— sind  coss  cos2^  -  ~sin2d  sins  sin20 

& 

Ra,  =cosd  coss  sin^  +  cos2d  sins  cos0 


For  an  explosive  source,  only  the  Rayleigh  wave  has  a  source  function: 


‘‘S  =0. 


*kur 


II.  Scattering  angular  integration: 

The  angular  integrations  for  determining  the  scattering  coefficient  (Equation  6) 
are 

©4  =  CoC2 

©a  “Ao»2  +  Coe?  +  2C0c2c0  +  Eo(aiCi+aoCa)  +  FgbgCj 
©o  =»Aoao  ■+•  Bobo  +  QjCo  +  Doaobo  +  EoSoO)  +  FoboCo 

where 


A0  —(vpw2)2  Bo=((2+v)X)2  Co  =((2+v)/x)2 

D0  =-2v(2+v)pXw2  Eo  =-2v(2+v)p/iw2  F0  =2(2+v)2Xp 

And  a’s,  b’s,  and  c’s  for  RR  conversion  are 

ao  «UIBU„ 

bo  -(k-unB-aiui01)  (k.u^-a.u,,,,) 

Ca  =2kmUtm  k„U„  Cl  =•( knU.n+a.U™,) ( kBU,B+d,Ura) 

Co  — 23,U1Bd,U,B 

for  RL, 

a.  -U^U,. 
b0  =0 

Ca  —  kBU  nnkgUfn 

for  LR, 

a.  =U#fflU„ 

bo-0 

ca  —kBU*BkBU„ 
and  for  LL, 

*i  ao  — 0 


ao  —0 

Ci  Cg  =0 


ao  =0 

Cl  —  dIUrBaiU#B+kBUimdIU  ,B  c0  =0 


c2  ^kmUfnjkaUin  Cj  — Co  — 0 

Appendix  II 

Locked  Mode  Approximation 

The  locked  mode  approximation  proposed  by  Harvey  (1981)  is  used  to  synthesize 
body  waves  by  surface  wave  normal  mode  theory.  A  rigid  cap  layer  is  placed  at  great 
depth  of  the  original  model.  This  acts  as  to  trap  the  leaking  modes  and  force  them  to 
become  the  locked  modes  (Watson,  1972).  From  normal  mode  theory,  the  imposition 
of  this  cap  layer  makes  the  branch  line  shrink  and  open  a  space  to  accept  the  leaking 
modes  as  they  migrate  from  lower  Riemann  surfaces  into  the  real  wavenumber  axis. 
Equivalently,  this  breaks  down  the  continuous  wavenumber  response  of  body  wave 
into  discrete  wavenumber  spectra  of  surface  waves. 

The  locked  mode  approximation  method  inherits  most  of  the  merits  of  surface- 
wave  normal  mode  theory.  The  advantages  of  this  method  include: 

(1)  The  layer  response  is  independent  of  the  source  or  receiver  positions.  A  data  file 
of  layer  response  can  be  set  up  first,  which  takes  most  of  computation  time.  The  data 
file  can  then  be  applied  to  many  different  source  and  receiver  configurations  without 
repeating  the  calculation  of  layer  responses. 

(2)  The  parameters  used  in  calculation  are  mostly  real  numbers  which  help  in  saving 
space  and  computation  time.  The  attenuation  effect  can  be  included  easily  by  using 
perturbation  theor  (Anderson  et.  al.  1965). 

(3)  The  method  can  handle  a  wide  range  of  frequencies,  say  0.001  to  1000  Hz,  as  well 
as  an  arbitrary  number  of  layers. 

(4)  All  of  the  rays  possibly  excited  in  the  model  are  automatically  included.  This 
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makes  the  method  appropriate  to  use  especially  for  complicated  models  and  large  dis¬ 
tances.  This  is  also  the  goal  generalized  ray  theory  attempts  to  reach. 

(5)  The  position  of  pole  (locked  mode)  is  exactly  known,  and  its  contribution  can  be 
correctly  calculated.  This  is  significant  compared  to  other  wavenumber  integration 
methods.  Those  methods  usually  use  complex  frequencies  to  shift  the  poles  away  from 
the  integration  path,  i.e.,  real  wavenumber  axis  before  taking  the  direct  integration 
(Bouchon,  1979).  Since  the  integrand  (layer  response)  varies  very  irregularly  along 
the  wavenumber  axis,  it  is  difficult,  if  not  impossible,  to  find  a  scheme  which  is  really 
efficient  to  sample  the  integration  path  especially  for  large  distances  or  high  frequen¬ 
cies.  The  locked  mode  theory  dose  not  have  this  problem.  Since  it  knows  where  the 
pole  is,  it  can  estimate  the  pole  contribution  exactly. 

To  make  the  ’locked  mode’  approximation  feasible,  several  requirement  need  be 
fulfilled  beforehand. 

(1)  All  of  the  dispersion  values  in  the  frequency  range  of  interest  must  be  found  com¬ 
pletely.  Insufficiently  determined  dispersion  values  always  produce  unsatisfactory 
results.  However,  the  number  of  locked  modes  generated  is  always  very  large,  say 
1000  modes.  Efficient  handling  of  this  large  number  of  modes  is  an  important  thing  in 
using  ’locked  mode’  approximation.  Wang  and  Herrmann  (1980)  proposed  a  method 
called  jumping  method  which  attempts  to  follow  each  pole  at  succeeding  frequencies 
rather  than  to  find  them  independently.  This  method  has  been  proved  successful 
especially  for  high  frequency  data.  Recently,  this  has  been  further  improved  to 
prevent  mode  missing  when  the  phase  velocity  of  a  mode  becomes  equal  to  a  layer 
velocity.  The  dispersion  curves  bend  abruptly  at  these  places  (Kerry,  1981).  But  this 
method  is  not  perfect.  The  modes  are  sometimes  jumped  over  in  certain  cases.  After 
some  effort,  we  finally  realize  that  there  may  not  exist  a  perfect  pole  searching  scheme 
which  can  fit  any  kind  of  problem  within  limited  computation  time.  To  solve  this 


problem,  an  independent  procedure  was  designed  to  ’patch  up’  the  missing  part  based 
on  the  dispersion  data  already  found. 

(2)  The  pole  position  (dispersion  value)  must  be  located  accurately.  In  a  recent  paper, 
Schwab  et  al.  (1984)  showed  that  dispersion  values  can  be  calculated  at  almost  unlim¬ 
ited  frequency  by  using  the  Knopoff  (1964)  matrix  method.  He  has  successfully 
tested  Harkrider’s  model  (Harkrider,  1970)  at  a  frequency  as  high  as  10,000  Hz.  We 
do  not  use  Knopoff’s  theory,  instead,  the  conventional  Dunkin  (1965)  compound 
matrix  is  applied  but  with  extensive  modification.  Using  the  same  model,  we  calculate 
the  dispersion  values  and  the  results  are  as  good  as  Schwab’s.  The  problem  of  accu¬ 
racy  of  dispersion  value  is  now  well  solved. 

(3)  The  terms  which  enter  the  final  solution  of  normal  mode  theory  need  also  be  cal¬ 
culated  accurately.  These  terms  include  eigenfunctions,  amplitude  factor,  ellipticity, 
and  attenuation  factor  (gamma).  Using  some  analytic  expressions  developed  by  Wang 
(1981)  following  Harkrider  (1980),  these  parameters  can  now  be  calculated  very  pre¬ 
cisely.  The  amplitude  factor  which  corresponds  to  the  pole  residue  is  a  sensitive 
parameter  difficult  to  determine.  Rather  than  taking  direct  derivative  (Harvey,  1981), 
it  is  calculated  by  the  energy  integrals  (Aki  and  Richards,  1980,  chap  7).  The  value  of 
the  Lagrangian  can  be  used  to  judge  the  accuracy  attained.  Our  results  up  to  now 
have  all  passed  this  test. 

There  seems  no  theoretical  difficulty  for  locked  mode  approximation  method 
now.  However,  this  method  has  an  internal  weakness.  That  is  where  to  put  the  cap 
layer  and  how  to  remove  the  unwanted  effects  caused  by  it.  Two  erroreous  phases 
might  be  generated  by  this  artificial  layer;  one  is  the  cap  layer  phase  which  is  caused  by 
real  arrivals  reflected  from  the  cap  layer,  and  the  other  is  the  truncation  phase  men¬ 
tioned  by  Harvey  (1981).  The  truncation  phase  arises  from  a  phase  velocity  filter 
which  is  used  to  suppress  the  cap  layer  phase.  To  avoid  a  cap  layer  phase  coming  too 
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early,  a  deeper  layer  should  be  used  for  larger  receiver  distance.  However,  a  deeper 
cap  layer  traps  more  waves  or  equivalently  creates  more  modes.  For  example,  455 
modes  for  200  km  cap  layer  becomes  1248  if  the  cap  layer  is  placed  at  800  km.  A 
compromise  must  be  set  to  conserve  computational  effort.  The  velocities  of  cap  layer 
must  also  be  defined  carefully.  They  need  be  large  to  yield  the  desired  body  waves  in 
the  crust,  but  not  too  large  to  generate  too  many  additional  poles.  Here  is  our  experi¬ 
ence,  which  can  just  be  used  as  a  rule  of  thumb.  If  Vc  is  the  velocity  of  cap  layer 
which  is  at  L  km  depth  and  V]  is  the  layer  velocity  above  the  cap  layer,  a  relation  can 
be  derived  assuming  simple  reflection: 


Ve 

L 


1+r 

1-r 

j/F 

l-r 


V, 

R 


r  represents  the  reflection  coefficient  desired  from  the  cap  layer  and  R  is  largest  dis¬ 
tance  to  be  used.  The  value  of  r  around  0.4  is  suggested,  which  makes  Vc  2.5  times 
lager  than  V,.  and  the  cap  layer  depth  grossly  equals  to  the  distance.  However,  if  the 
phase  velocity  filtering  is  done  properly,  the  distance  range  can  be  extended  (at  least 
doubled). 


A  proper  phase  velocity  filter  can  partly  smooth  out  the  unwanted  cap  layer 
phase,  but  not  entirely.  However  the  phase  velocity  filter  itself  also  itroduces  a  noise 
caused  by  the  trunction  of  filter.  The  truncation  phase  may  change  position  and  shape 
in  the  synthetic  seismograms  using  different  phase  velocity  filter  values.  Hence,  it  can 
be  identified  by  comparing  the  seismograms  obtained  from  different  phase  velocity 
filters.  A  frequency  varying  phase  velocity  filter,  sometimes,  is  found  useful  in 
suppressing  the  truncation  phase. 

Figure  22  and  23  show  some  synthetic  seismograms  made  from  the  locked  mode 
approximation.  In  Figure  22,  the  tangential  component  at  distances  from  100  too  500 
km  is  generated  using  a  simplified  one  layer  over  halfspace  model  (SCM,  Wang  and 
Herrmann,  1980).  The  arrival  time  of  each  phase  can  be  exactly  identified  Figure  23 


-  32  - 


displays  ten  basic  seismograms  for  the  CUS  model  from  double-couple  and  explosive 
source  at  a  distance  of  200  km.  These  seismogram  are  really  of  high  quality. 


-  34  - 


G.S. 

D.E. 

(m) 

K) 

TABLE  2 

Q  Values  of  Synthetic  Coda 
(Theoretical  Q  =  1000) 


DIST=100Km 


DIST=300Km 


5 


G.S.  **  Geometric  Spreading  Factor  »  l/r“ 

O.E.  Duration  Expansion  =*  a^r 

DIST  =“  Scattering  Distance 

dtp  —  time  interval  for  succeeding  pulses. 

In  each  pair,  the  upper  number  is  for  coda  between 
0  and  100  second  and  the  lower  one  is  for  100  and 
300  second. 
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TABLE  3 

Q  Values  of  Synthetic  Coda 
(Theoretical  Q  =  150) 


DIST=100Km 


DIST=300Km 


5 

10 

dt 

p-. 
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Figure  1  Observed  seismograms  from  ECTN  network.  (a)  is  the  original  trace,  (b)  is 
the  rectified  trace,  (c)  is  the  trace  narrow  bandpass  filtered  at  1  Hz.  The  coda  seems 
being  composed  of  many  small  ’pulses’  The  coda  envelope  is  not  so  uniform  if 
examining  closely.  The  horizontal  bar  indicates  a  60  second  time  window. 


Figure  3.  A  comparison  of  waveforms  of  direct  surface  waves  and  the  scattered  waves. 
The  numbers  at  the  end  of  each  seismogram  indicates  the  source-receiver  distance  for 
surface  waves  and  the  source-  scatterer-receiver  distance  for  the  scattered  waves. 


110  km 


300  km 


Figure  4.  Wave  type  conversion  test.  The  upper  trace  in  each  pair  of  seismograms  does 
not  include  wave  type  conversion  and  the  lower  one  does.  The  tests  are  for  scatterers 
at  the  distance  110  km  and  300  km,  and  at  the  angles  0,  90  ’  150  degrees  on  an 

ellipse. 


Figure  5.  Mode  conversion  test.  Cases  discussed  are  for  mode  conversion  permitted 
up  to  1,  5  and  10  modes.  Three  scattering  distances  110  km,  300  km  and  500  km  are 
included.  Two  columns  correspond  to  cases  with  or  without  wave  type  conversion 
respectively.  Plote  a,  b  and  c  are  for  three  different  components. 
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Ficrure  8.  Mode  conversion  test  for  an  explosive  source.  The  scattering  distance  is  at 
110  km.  At  this  distance,  both  mode  conversion  and  wave  type  conversion  are 

important. 
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Figure  7.  The  effect  of  scattering  at  different  location  of  an  ellipse.  Three  traces  at 
each  location  are  for  the  scatterers  at  the  distances  110  km,  300  km  and  500  km  with 
the  source-receiver  distance  at  100  km.  The  amplitude  of  each  trace  has  been 
normalized.  The  time  span  covered  is  30  seconds.  Plots  a,  b  and  c  are  for  three 
different  components. 
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Figure  8.  The  effect  of  scatterer  depth  fa)  shows  .  i 

km  scattering  distance  and  (b)  is  vertical  com  non  Jn?  ta“,gentlaI  component  with  200 
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Figure  9.  The  effect  of  different  dej 
distance  is  110  km  and  the  scatterer 
distance  is  100  km.  High  scattering 
components  and  broadens  the  scatterin 
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■ee  of  scattering  attenuation.  The  scattering 
s  at  the  depth  of  5  km.  The  source-receiver 
attenuation  partially  absorbs  high  frequency 
pulse. 
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Figure  11a.  The  scattering  waveforms  and  their  spectra  under  different  attenuation 
conditions.  The  scattering  distance  is  200  km  and  the  scatterer  is  5  km  deep  in  CUS 
model.  The  high  frequency  components  are  attenuated  equally  in  both  cases.  The 
spectra  plotted  are  in  the  same  sequence  as  the  seismograms  but  shifted  by  a  factor  of 
ten  for  clarity.  Their  high  frequency  spectra  slopes  are  also  indicated. 
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Figure  lib.  Similar  as  for  Figure  11a,  but  for  scatterer  at  15  km  deep.  Scattering  Q 
absorbs  the  high  frequencies  in  the  main  pulse.  The  anelastic  attenuation  is  weaker 
because  of  the  deeper  scatterer. 
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Figure  12a.  A  comparable  case  to  Figure  11.  The  scattering  distance  is  500  km  now. 
After  traveling  large  distance,  the  high  frequencies  are  absorbed  especially  for  anelastic 
attenuation.  In  this  case,  the  anelastic  Q  controls  the  decay  of  wave  energy. 
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Figure  12b.  Same  as  for  Figure  12a,  but  for  scatterer  at  15  km  deep.  Note  that  the 
high  frequency  spectrum  slope  difference  between  OUOQ  ( undissipative)  and  10UOQ 
(scattering  attenuation)  cases  are  very  similar  for  all  cases  in  Figures  11  and  12.  The 
scattering  Q  is  dominant  in  controlling  the  frequency  decay. 
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‘  i  Figures  13.  The  geometrical  spreading  factor  for  scattering  waves  in  the  layered 

medium.  Displays  are  for  Z  component  (a)  and  SH  component  (b)  with  epicentral 
[  distance  at  100  km  in  CUS  model.  The  scattering  distance  varies  from  110  km  to  600 

km.  The  top  curve  is  from  the  average  peak-to-peak  amplitude  of  composed  scattering 
^  waveforms,  and  the  lower  two  curves  are  for  amplitudes  after  narrow  bandpass  filtered 

r  at  1  Hz  and  3  Hz  respectively.  No  Q  has  been  applied.  Because  of  mode  conversion, 

[■-  we  can  see  much  fast  amplitude  decay  in  the  front  part  (  <200km). 
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Figure  14.  Scattering  geometrical  spreading  in  the  WUS  model.  Two  components, 
vertical  and  tangential,  are  displayed,  (a)  is  for  source  -receiver  distance  at  300  km 
and  (b)  at  100  km.  Note  that  the  slope  changes  at  600  km  for  (a)  and  at  200  km  for 
(b).  The  WUS  model  has  more  rapid  geometrical  decay  rate  than  the  CUS  model. 
The  slope  difference  between  the  two  regions  is  also  larger. 
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Figure  15.  A  modification  of  geometrical  spreading  correction  of  Sato’s  surface-wave 
scattering  model.  A  higher  decay  rate  of  t"06  to  t_l  due  to  mode  conversion  and  wave- 
type  conversion  is  assumed  for  fhe  lapse  time  ratio  smaller  than  2. 
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Figure  16b.  The  duration  change  measured  from  Figure  16a  are  plotted  in  linear  scale 
with  the  distance.  A  duration  expansion  rate  0.05  with  distance  is  obtained.  The 
upper  and  lower  sets  of  data  are  from  SH  and  Z  synthetics,  respectively. 
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Figure  17a.  Superposition  of  scattering  pulses  to  form  the  coda  envelope.  This 
simulates  the  eastern  United  States  case.  Spreading  factor  m  largely  affects  the  coda 
decay.  Also  note  that  the  analysis  of  first  part  (I)  and  the  rest  (II)  yields  different 
apparent  Q  values. 


Figure  17b.  Same  as  for  Figure  17a.  This  simulates  the  western  United  States  case. 
Since  the  pulse  superposition  has  been  saturated  by  the  rapid  decay  due  to  Q,  the 
inferred  Q  is  not  dependent  upon  the  shape  factors  of  the  individual  scattered  pulses 

much. 


Figure  18.  A  comparison  of  body  wave  and  surface  wave  scattering.  The  scattering 
distance  is  110  km  with  epicentral  distance  at  100  km  in  the  CUS  model.  Some  body 
and  surface  wave  scattering  interaction  can  be  seen. 
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Figure  19.  Body-wave  scattering  expanded  at  three  scattering  distances.  The  scattered 
waveform  is  complicated  and  contains  high  frequencies.  P-wave  scattering  is  strong  on 
the  Z  component,  and  S-wave  scattering  is  strong  on  the  two  horizontal  components. 
The  complicated  waveforms  show  the  ray  path  effect  in  the  layered  structure. 


Figure  20a.  The  spectra  of  body-wave  and  surface-wave  scattered  pulses.  Three  traces 
are  for  complete,  body-wave  only,  and  surface-wave  only.  The  spectra  of  surface- 
wave  only  is  shifted  down  by  a  factor  of  100  for  clarity.  Note  that  the  spectra  of  the 
scattered  body  waves  is  enriched  in  high  frequencies,  even  though  no  attenuation  is 
included. 
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Figure  20b.  Same  as  for  Figure  20a  with  a  scattering  Q  included.  The  attenuation  of 
body-wave  scattered  waveforms  is  less  than  the  surface- wave  scattered  arrivals. 
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The  ,decay  of  hl®hL  fre^uency  component  under  different  attenuation 
conditions  for  surface-wave  and  body-wave  scattering.  The  body-wave  scattering  shows 
attenuatton^ compared  to  surface-wave  scattering  attenuation.  Each 
spectrum  has  been  shifted  by  a  factor  of  ten  for  clarity. 
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Figure  21b.  same  as  for  Figure  21a,  but  for  a  scattering  distance  of  500  km.  The 
frequency  content  of  body-wave  scattering  seems  not  affected  as  much  by  the 
attenuation. 


Figure  22.  The  synthetic  seismograms  from  ’locked  mode’  approximation.  SH 
component  from  a  single  layer  over  hrdfspace  model  (SCM)  is  plotted  with  different 
epicentral  distance.  The  arrivals  from  aifferent  reverberation  are  followed  through  the 
section. 
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Figure  23.  Ten  basic  type  complete  seismograms  generated  by  locked  mode 
approximation  method.  The  seismograms  rather  have  excellent  quality. 
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